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Abstract 


The common Asian toad Duttaphrynus melanostictus (Schneider, 1799) complex has a wide distribution ranging from western foot- 
hills of the Himalaya to the easternmost range of the Wallacea, with the evidence of human-mediated introductions to some other ar- 
eas. In the entire distribution range, the complex is formed by several evolutionary clades, distributed mostly in South-East Asia with 
unresolved taxonomy. In the northwestern edge of its distribution (Pakistan), the name D. melanostictus hazarensis (Khan, 2001) has 
been assigned to local populations but its biological basis remained, so far, understudied and unvalidated. Therefore, we re-evaluated 
the available genetic data (mitochondrial and nuclear) to show the relationships between Pakistani populations (including the type 
locality of D. m. hazarensis) and others from across the range. Our results showed that Pakistani populations are associated with 
one, deeply diverged, well-supported and widely distributed clade (so-called Duttaphrynus sp. 1 according to /6S, or clade B based 
on tRNAGIy-ND3), that has already been detected in previous studies. This clade is further distributed in India, Nepal, Bangladesh, 
Malaysia, Singapore, and Indonesia and is characterized by a low level of genetic variability. This further suggests that both natural, 
as well as potential human-mediated dispersal, might have played an important role in setting up the current phylogeographic and 
distribution pattern of this clade. The clade is deeply divergent from other clades of the complex and represents a taxonomically 
unresolved entity. We here argue that the clade Duttaphrynus sp. 1/B represents a distinct species for which the name Duttaphrynus 
bengalensis (Daudin, 1802) comb. nov. is applicable, while the description of D. m. hazarensis does not satisfy the rules of the In- 
ternational Code of Zoological Nomenclature. 
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Introduction 


The distribution of the Duttaphrynus melanostictus 
(Schneider, 1799) complex (Dubois and Ohler 1999) 
more or less corresponds to the area of the Oriental 
zoogeographic realm. It is a widespread group of toads 
with a geographic range throughout South and South- 
East Asia, southern China, including islands of Taiwan, 


Hainan, Sumatra, Java, Borneo, and Sumba (Frost 
2021). The group has a high potential to spread into 
new localities and was artificially introduced to many 
regions of Asia as well as out of Asia (see Reilly et al. 
2017; Vences et al. 2017; Othman et al. 2020; Soorae et 
al. 2020). These introductions are well-monitored be- 
cause they can pose a Serious threat to the local native 
fauna (Pettit et al. 2021). 
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The Indian subcontinent and northern Pakistan form 
the western border of the D. melanostictus complex 
natural range where its distribution is limited to the 
edge of Shivalik Hills, i.e. to several districts of Khyber 
Pakhtunkhwa, Punjab, and Azad Jammu and Kashmir 
Provinces (Mertens 1969; Khan 1972, 2001; Masroor 
2012). Also, the presence of the complex in Sindh 
Province cannot be excluded (Masroor 2012) as this 
toad is known from several localities in Gujarat (India) 
close to the Pakistani border (Sarkar 1984). Such a 
distribution pattern signifies that D. melanostictus has 
been recorded in both the Palearctic and Oriental region, 
which are the main zoogeographic realms of Pakistan, 
(see Masroor 2012: fig. 20). Surprisingly, populations in 
the northwestern part of the range have been described 
as. Duttaphrynus melanostictus hazarensis (Khan, 
2001) based on morphological data. According to Khan 
(2001), specimens of these populations differ from other 
Asian populations of D. melanostictus by bearing bony 
ridge, kidney-shaped parotids, subarticular tubercles 
that are doubled under the penultimate phalanx of all 
fingers, the absence of rostral ridge, and a light brown 
dorsal coloration. However, the taxon description and 
comparison provided by Khan (2001) are very brief. 
Further, the original description neither provides a 
voucher nor a catalogue number for the holotype, nor 
information of the type collection or museum designation. 
Therefore, it does not fulfil the International Code of 
Zoological Nomenclature (Art. 16.4.2, ICZN), rendering 
the subspecific name D. m. hazarensis (Khan, 2001) 
invalid. In order to simplify the terminology in our text, 
we will, however, refer to the Pakistani population and 
related lineage as “hazarensis’’. 

Due to the topographic heterogeneity and habitat 
diversity, it might be plausible that Pakistani amphibian 
populations of different genera with the Oriental affiliation 
may hold a different evolutionary history compared to 
their populations from other parts of South-East Asia. The 
significance of the far-western part of the Himalaya and the 
Indo-Gangetic plain for understanding the biogeography 
of taxa distributed across these geomorphological 
domains has been recently shown in several species 
groups (see Allopaa, Chrysopaa, Microhyla, Nanorana, 
Scutiger, Sphaerotheca;, Hofmann et al. 2017, 2019, 2021; 
Jablonski et al. 2020, 2021). However, except for Khan’s 
(2001) original work, populations of D. melanostictus 
from Pakistan, including toads from the type locality 
where the “hazarensis” population has been described, 
have not been studied so far in a broader biogeographical 
context (cf. Hussain et al. 2020; Akram et al. 2021). 
According to Khan (2001), the type locality has a wider 
range and is defined as “Ooghi, Mansehra and Datta, 
Hazara Division, eastern NWFP, Pakistan’, 1.e., covering 
several villages and districts in northern Pakistan where 
these toads occur. Recent molecular studies focusing on 
the genus Duttaphrynus were based on dense sampling, 
covering a wide range of South-East Asia (mostly), and 
different DNA markers for D. me/anostictus (see Hasan et 
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al. 2014; Wogan et al. 2016; Reilly et al. 2017; Vences et 
al. 2017; Othman et al. 2020; Soorae et al. 2020; Bisht et 
al. 2021). Because the resulting DNA sequence variation 
showed deep genetic structuring, we reviewed and 
reanalysed these data together with new sequence data 
of D. melanostictus from different localities of Pakistan, 
including the locality of “hazarensis”’ (Datta, Khyber 
Pakhtunkhwa; 34.2962°N, 73.2649°E). Our main goals 
of the study are to (1) find out the genetic affiliation of the 
Pakistani population of D. melanostictus in the context 
of known molecular phylogeny of the complex, (i1) show 
phylogeographic relationships of the western populations 
to others inside the detected clade, and (iil) comment 
on the validity and taxonomy of the Pakistani taxon in 
light of molecular phylogeny and nomenclature of the 
D. melanostictus complex. 


Materials and methods 
Sampling and DNA extraction 


Tissue samples of six individuals of the Duttaphrynus 
melanostictus complex from northern Pakistan were col- 
lected from five different localities of the Khyber Pakh- 
tunkhwa Province during the field studies between 2018 
and 2019. Individuals were documented photographically 
by DJ (Nikon D810, and 105 mm macro lens) and sam- 
pled for toe or blood before release or being deposited in 
scientific collections (Suppl. material 1: Table S1). Tissue 
samples were transferred into 96% ethanol, kept at ambi- 
ent temperature during the time of the fieldwork, and later 
stored at -20 °C. DNA was isolated using the DNeasy Tis- 
sue Kit (QIAGEN), following the manufacturer’s proto- 
cols, and stored at -20 °C. 


Amplification and sequencing 


We amplified a fragment of the nuclear (nDNA) rhodop- 
sin (Rho) gene (~316 bp), of the /6S rRNA (~560 bp), and 
approximately 469 bp of mitochondrial DNA (mtDNA) 
that encodes part of the cytochrome oxidase c subunit III 
(COITD, tRNA glycine (tRNAGIy), and part of the NADH 
dehydrogenase subunit 3 (VD3) using the primers and 
PCR conditions following previous studies (Vences et al. 
2017; Hofmann et al. 2019; Bisht et al. 2021). PCR prod- 
ucts were Sanger-sequenced in both directions by Macro- 
gen Europe (Amsterdam, The Netherlands; http://www. 
macrogen-europe.com) and sequence data edited where 
necessary in Sequencher 5.4 by eye. 


Dataset and phylogenetic analyses 


We complemented our new mitochondrial and nuclear 
DNA sequences dataset with published and unpublished 
sequence data from GenBank (for details see Suppl. 
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material 1: Table S1). We mostly followed datasets pro- 
vided by Wogan et al. (2016), Vences et al. (2017), and 
Bisht et al. (2021) as well as others with available se- 
quences for D. melanostictus. To map the phylogeograph- 
ic pattern of detected clades, we used published georef- 
erenced localities of the sequence origin. If coordinates 
were not available for sequences, we georeferenced them 
according to the provided name of the locality or speci- 
men voucher numbers connected to a particular museum 
collection. In some cases, sequences could not be georef- 
erenced (e.g., “India”; especially for /6S dataset) or with 
uncertainty only (e.g., “Uttarakhand”, “Kalimantan”; de- 
tails in Suppl. material 1: Table S1). For an overall view 
on the distribution of the complex we downloaded data 
for D. melanostictus (~20,000 records) available in Glob- 
al Biodiversity Information Facility (GBIF, accessible via 
www.gbif.org). 

New sequences were aligned in MEGA X v.11.0 
(Kumar et al. 2018). We followed Vences et al. (2017) 
and trimmed the mitochondrial sequences, retaining 411 
bp corresponding to 71 bp of tRNA-Gly and 340 bp of 
ND3. We excluded the following three /6S sequences of 
the data set of Bisht et al. (2021) due to low data quali- 
ty and potential sequencing errors: EU071741, U52798, 
EU071765. Also, one sequence (AY458592) was ex- 
cluded from the data set of Vences et al. (2017) due to 
uncertain geographic origin. The resulting data sets con- 
sisted of 150 sequences for /6S (560 bp), 496 sequences 
for tRNAGIy-ND3 (411 bp), and ten sequences for Rho 
(300 bp) (Suppl. material 1: Table S1). The latter sequence 
alignment was used only for the haplotype network con- 
struction. All new sequences were deposited in GenBank 
(see Suppl. material 1: Table S1 for accession numbers). 

Bayesian phylogenetic analyses were carried out 
with MrBayes v.3.2.6 (Ronquist et al. 2012) using the 
best-fit model of sequence evolution as determined for 
16S (GTR+I+G) by jModeltest v.2.1.10 (Darriba et al. 
2012) and tRNAGI/y-ND3 (four partitions: K80, K80+G, 
HKY+G, GTR+I+G) by PartitionFinder v.1.1.1 (Lanfear 
et al. 2012). We run MrBayes with a random starting 
tree for five million generations using four parallel Mar- 
kov chain Monte Carlo (MCMC) simulations with four 
chains, sampling trees every 500 generations, until reach- 
ing an average standard deviation of split frequencies 
of <0.01. Stationarity and convergence of the runs were 
confirmed using the software Tracer v.1.7.2 (Rambaut et 
al. 2018). The first 25% of the samples of each run was 
discarded as burn-in. 

For reasons of comparison, we also inferred the BI tree 
for /6S based on the secondary structure of the ribosom- 
al RNA using RNAsalsa 0.8.1 (Stocsits et al. 2009). As 
the initial input, we used the constraint file based on the 
secondary structures of Bos taurus 16S rRNA which 1s 
provided with the RNAsalsa package. The sequence was 
then aligned with our /6S data using the MUSCLE algo- 
rithm in MEGA X. Starting with this initial alignment and 
the constraint file, RNAsalsa predicts RNA secondary 
structure and enhances structural alignment that results 
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in a final multiple sequence alignment together with a 
consensus structure. Unambiguous stem pairs were then 
derived based on the consensus structure from RNAsalsa, 
specified in the MrBayes input file, and analysed under 
the doublet model (16x16) proposed by Schoninger and 
von Haeseler (1999). 

To assess the evolutionary distance between our sam- 
ples and the Duttaphrynus sister clades recovered from 
the phylogenetic analysis, we calculated the uncorrected 
p-distances for the /6S rRNA and tRNAGIy-ND3 data set 
separately using Mega-X based on 1,000 bootstrap repli- 
cations, with the pairwise deletion option, and by consid- 
ering both transitions and transversions. For clades that 
correspond with the Pakistani population of D. melanost- 
ictus, and for the Rho sequence data we created a medi- 
an-joining haplotype network using PopART (Leigh and 
Bryant 2015). 


Results 


The newly generated sequences of /6S, tRNAGIy-ND3, 
and Rho from Pakistan belong to one, well-support- 
ed, monophyletic group in all investigated markers 
(Figs 1-3). This clade is known from previous studies, 
particularly Wogan et al. (2016), Vences et al. (2017) and 
Bisht et al. (2021). 

For /6S, populations from Pakistan correspond with 
the clade called Duttaphrynus sp. 1 (sensu Bisht et al. 
[2021]) that includes populations from Bangladesh, 
India, Nepal, and surprisingly, Malaysia (continental 
as well as Borneo part; Figs 1A, 2A). This clade has a 
sister position to clade Duttaphrynus sp. 2 (sensu Bisht 
et al. 2021) including populations from China, Laos, 
Malaysia, Myanmar, and Vietnam. Overall, seven clades 
of the D. melanostictus complex were detected in the /6S 
phylogeny. Uncorrected p-distances between these clades 
range from 2.3 to 4.9% (Table 1). The genetic distance 
between Duttaphrynus sp. 2 and Duttaphrynus sp. 1 
clades is 2.7%. The clade Duttaphrynus sp. 1 1s formed by 
12 haplotypes in 24 sequences with low internal sequence 
variation (m=1.0%). Sequences from Pakistan have a 
closer relationship (one mutation step) to sequences from 
peninsular India (KU589218, EU071749, EU071752, 
EU071754, EU071757, KY000453-KY000456) and 
Nepal (MF319517; Fig. 1B). Other sequences (haplotypes) 
differ by two to six mutation steps. Surprisingly, two 
sequences from Malaysia (including Borneo; AB530648, 
AB331714) are separated by only two mutation steps 
from Pakistani sequences, or are even identical to 
some sequences from India or Nepal (Figs 1B, 2A). 
Significantly, some sequences from India could not be 
georeferenced and, thus, not displayed on the map in Fig. 
2A (see Suppl. material 1: Table S1). However, the overall 
phylogeographic pattern suggests a distribution of the 
clade Duttaphrynus sp. 1 in the northern part of the Indian 
subcontinent with a potential range limit in Myanmar, 
where it is replaced by the sister clade Dutaphrynus sp. 2. 
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Figure 1. The molecular phylogenetic reconstruction of the Duttaphrynus melanostictus complex on 16S rRNA (A) sensu terminol- 
ogy of Bisht et al. (2021) and tRNAGIy—-ND3 (C) sensu Vences et al. (2017) together with median-joining networks (B, D) calculated 
from clades that include sequences of the “hazarensis” toads from Pakistan, now D. bengalensis comb. nov. (Daudin, 1802). The 
numbers on branches represent Bayesian Posterior Probabilities/bootstraps showing branch support. Sequence details are presented 
in Suppl. material 1: Table $1. The same colors of clades used for particular mtDNA trees reflect clade affiliation. The pictured 
individual represents a population from Shahdara, Margalla Hills, Islamabad, Pakistan. 


Populations of this ‘sp. 1’ clade in Western Malaysia and 
Borneo (Sarawak) might have been artificially introduced 
(see Discussion below and Fig. 2A). 

For tRNAGIy-ND3, Pakistani sequences correspond to 
the clade B (sensu Vences et al. [2017]) that includes pop- 
ulations from Indonesia (Kalimantan), West Malaysia, 
and Singapore (Figs 1C, 2B). This clade has a sister rela- 
tionship with the clade C distributed in Myanmar. Over- 
all, we confirmed seven clades of the D. melanostictus 
complex in the phylogeny of tRNAG/y-ND3. Uncorrected 
genetic distances between clades range from 2.7 to 9.3% 
(Table 2). The genetic distance between closely related 
clade B and C is 5.6%. The clade B includes three haplo- 
types. Two of them are found in Pakistan and differenti- 
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ated by one mutation step from each other. The third hap- 
lotype differs by four mutation steps from the Pakistani 
sequences and can be found in Indonesia (KU183404), 
Malaysia (KU183466, KU183468-KU 183470), and Sin- 
gapore (KU183475, KU183476, KU183484). Internal 
sequence variation in clade B 1s very low (a = 0.3%). 
Insufficient sampling across the Indian subcontinent did 
not allow for assessing range limits of the clade B in rela- 
tion to other clades in the resulting phylogeny. However, 
detection of the clade B in Singapore, West Malaysia, and 
Indonesia, (Borneo; Kalimantan, Fig. 2B) suggests an ar- 
tificial, probably human-mediated introduction. 

The results on rhodopsin (Rho) confirm mitochon- 
drial data. It shows that D. melanostictus from Pakistan 
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Table 1. Uncorrected genetic distances (%) between investigated /6S clades of the Duttaphrynus melanostictus complex as shown 
in the Fig. 1A. Given are p-distance values (lower-left) and standard error estimates (upper-right). 
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Table 2. Uncorrected genetic distances (%) between investi- 
gated tRNAGIy-ND3 clades of the Duttaphrynus melanostictus 
complex as shown in the Fig. 1C. Given are p-distance values 
(lower-left) and standard error estimates (upper-right). 
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represents an independent group (Fig. 3). This group 
is formed by sequences from Pakistan and differs by 
three mutation steps from populations from Vietnam 
(DQ283967; based on the geographic position it proba- 
bly represents the mtDNA clade Duttaphrynus sp. 2 in 
16S, or clade G in tRNAGIy-ND3, respectively) and by 
two to six mutation steps from D. melanostictus from 
India (AF249097; without precise locality). The Indian 
sequence cannot be related to any known mitochondrial 
phylogenetic clade due to insufficient sampling and un- 
known geographic origin. The uncorrected genetic dis- 
tances between Rho sequences of the D. melanostictus 
complex range between 0.7 and 1.3%. 


Discussion 


The most prominent results of our study are the following: 
(i) the “hazarensis” population from Pakistan belongs to 
a widely distributed and taxonomically so far unresolved 
clade (see below), (11) Pakistani mitochondrial sequences 
represent unique haplotypes with a close relationship to 
other populations from northern parts of the Indian sub- 
continent, (111) comparing the overall phylogeographic 
pattern, the clade Duttaphrynus sp. 1/clade B has been 
most likely artificially translocated to South-East Asia, 
probably from the Indian subcontinent, (iv) the name for 
the clade Duttaphrynus sp. 1/clade B is available and, 
thus, should be applied (see Taxonomy). These results are 
thus important for both, evolution and taxonomy of the 
D. melanostictus complex. However, they still underline 


previous conclusions that the complex needs comprehen- 
sive biogeographical research with taxonomic revision 
based on dense DNA sampling, especially from the In- 
dian subcontinent (see Wogan et al. 2016; Vences et al. 
2017; Bisht et al. 2021). 


Genetic affiliation and origin of Pakistani 
populations 


By re-analysing /6S sequences of the D. melanostictus 
complex we showed a pattern of three main, geograph- 
ically widely distributed clades (called as Duttaphrynus 
melanostictus, Duttaphrynus sp. 1, and sp. 2 sensu Bisht 
et al. 2021) with potential distribution boundary between 
them in western Myanmar and Malay Peninsula/Sumatra 
(Fig. 2). Populations from Pakistan belong to the Dut- 
taphrynus sp. 1 clade that is distributed in the northern 
part of the Indian subcontinent (Fig. 1; Bisht et al. 2021). 
Unfortunately, limited sampling in central India does not 
allow delimitating the range of this clade in the south- 
ern parts of India. The Indian subcontinent 1s, however, 
absolutely crucial to understand the range of this clade 
as southern India is inhabited by different evolutionary 
clades of the complex, 1.e. D. brevirostris, D. microtym- 
panum, and D. melanostictus (sensu Bisht et al. 2021; 
Figs 1, 2). Re-evaluation of tRNAG/y-ND3 sequence data 
of the complex revealed similar results as for /6S, show- 
ing two widely distributed clades (A and G) and four geo- 
graphically disjunct clades assigned to the northern and 
western parts of Indochina (C, E, D, F; Wogan et al. 2016; 
Vences et al. 2017; Figs 1, 2). Clade B (sensu Vences et 
al. 2017 as a part of the “mainland clade” sensu Wogan et 
al. 2016), where Pakistani toads belong to, was previous- 
ly reported from Indonesia (Kalimantan, Borneo), Ma- 
laysia, and Singapore (Wogan et al. 2016; Vences et al. 
2017). However, since these former studies did not focus 
on this specific clade, no conclusions were drawn about 
its distribution area. Despite limited data from the Indian 
subcontinent, the genetic clustering of Pakistani samples 
within clade B, and the similar phylogenetic placement 
observed in the /6S tree suggest a wide distribution range 
of this clade. The low level of intra-clade sequence differ- 
entiation and star-like pattern of the haplotype networks 
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Figure 2. Phylogeographic reconstruction and distribution of major clades of the Duttaphrynus melanostictus complex across the 
species range based on /6S rRNA (A) and tRNAG/Iy—ND3 (B) sequences available in GenBank and from own data (see Suppl. ma- 
terial 1: Table S1). The clade colors follow those of Fig. 1, clade terminology is derived from Bisht et al. (2021; 76S rRNA), and 
Vences et al. (2017; tRNAGIy—-ND3). The pictured individual represents the population of Duttaphrynus bengalensis comb. nov. 
(Daudin, 1802) from Datta, Mansehra, Pakistan. Coloration reflect the main mtDNA clades according to Fig. 1. Distribution records 
(white circles) have been downloaded from GBIF (gbif.org). The red line (panel A) depicts the area considered as “Bengal”, the 
type locality of Bufo bengalensis. 


implies a very close relationship between geographically For the clade Duttaphrynus sp. 1/clade B, we revealed a 
distant populations and potentially rapid expansions to _phylogeographic structure resembling patterns detected in 
wider areas, supplemented by recent human-mediated in- other widely distributed amphibians of the Indo-Gangetic 
troductions (see below). plain. These species colonized Pakistan from the Oriental 
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Figure 3. Median-joining alleles network of sequences from the Duttaphrynus melanostictus complex and closely related members 
of western Palearctic Bufonidae based on rhodopsin (Rho) gene available in GenBank and from own data. The colors correspond to 
those used in Fig. 1. Note, the origin of sequence AF249097 from India cannot be assigned due to unknown geographic details. The 
pictured individual represents a population of Duttaphrynus bengalensis comb. nov. (Daudin, 1802) from Datta, Mansehra, Pakistan. 


region, namely Microhyla nilphamariensis Howlader, Nair, 
Gopalan, Merila, 2015 and Sphaerotheca maskeyi (Schle- 
ich & Anders, 1998) (Jablonski et al. 2020, 2021). Both 
Species are wide-ranging, occupying central and north- 
ern parts of the Indian subcontinent, from eastern India, 
through Bangladesh, Nepal, central, northern, and west- 
ern India to Pakistan. Their distribution more or less fol- 
lows the Indo-Gangetic plain including the Brahmaputra, 
Yamuna, Ganga, and Indus rivers. Thus, we assume that 
the clade Duttaphrynus sp. 1/clade B colonized new areas 
along these riverine corridors (probably from the northern 
India and Nepal area; see the structure of the haplotype net- 
work Fig. 1B) and that this clade can also be found tn cen- 
tral and western India, eventually also southern Pakistan 
(see Sarkar 1984). However, as shown by the model of the 
potential distribution of the complex in Asia (see supple- 
mentary fig. S2 in Vences et al. 2017), the distribution of 
these toads might be scattered in central India and denser 
to the region of the Shivalik Hills (see also the georefer- 
enced records from the GBIF database Fig. 2). The Shiva- 
liks (also known as the Churia Hills) may provide optimal 
environmental conditions for these toads (or this specific 
evolutionary clade). It is a mountain range of the outer Hi- 
malayas that stretches from Pakistan about 2,400 km east- 
wards (10-50 km wide) close to the Brahmaputra River, 
spanning across the northern parts of the Indian subconti- 
nent. As the western edge of these hills is currently inhabit- 
ed by the clade Duttaphrynus sp. 1/clade B we can assume 
that such an environmental corridor (the Siwalik Corridor 
sensu Khan 2006) also affects the distribution of this spe- 
cies complex, enabling migration from the east to the west. 

Recently, two basic biogeographic scenarios have been 
hypothesized for the origin of different amphibian taxa 
distributed in northern Pakistan: 1) dispersal of ancestral 


lineages during the Miocene from East or Southeast China 
into the far northwestern part of the Himalaya-Tibet orogen 
along a warm temperate corridor that existed in significant 
parts of Paleo-Tibet’s interior (Hofmann et al. 2021), this 
corridor played probably an important role in the initial 
diversification process of Himalayan faunal organisms 
(Schmidt et al. 2012), (11) dispersion on widely distributed 
species characterized by the low level of genetic divergence 
that followed Shivalik Hills and the Indo-Gangetic plain 
and reached the territory of Pakistan (Jablonski et al. 
2020, 2021) as a result of the Pleistocene or Holocene 
colonisations (see Othman et al. 2020). Although these are 
presumable biogeographic models, the second scenario 
could be now applied for D. melanostictus from Pakistan. 
It is well known that the D. melanostictus complex is 
considered an invasive organism, with several evolutionary 
clades that were introduced out of its native range (Reilly 
et al. 2017; Vences et al. 2017; Othman et al. 2020; Soorae 
et al. 2020; Pettit et al. 2021). Evidence for potentially 
human-mediated introductions is visible from our anal- 
ysed data, showing large geographic distances between 
the source and apparently introduced populations (Fig. 2). 
Consistent with previous studies (Reilly et al. 2017; Vences 
et al. 2017; Othman et al. 2020; Soorae et al. 2020), proba- 
bly clades A and G were sources of introduced populations 
(sensu Vences et al. 2017, Fig. 1). These clades are widely 
distributed and covered by dense sampling. Moreover, these 
clades are distributed also in and around populated cities 
and ports of South-East Asia from where the artificial, hu- 
man-meditated introduction to other areas is highly evident 
(even to Madagascar or the United Arab Emirates; Vences 
et al. 2017; Soorae et al. 2020). On the other hand, limited 
sampling of the clade B (sensu Vences et al. 2017 which is 
equal to clade Duttaphrynus sp. 1 sensu Bisht et al. [2021] 
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in the /6S dataset) did not allow to evaluate the origin of 
this clade in the Malay Peninsula or Borneo. Based on the 
combination of /6S and tRNAGIy-ND3 data and the dis- 
tribution pattern of both detected clades (Duttaphrynus sp. 
1/clade B) we can hypothesize that Malaysian, Singapore, 
and Indonesian populations are results of the expansion of 
the range during ice-ages followed by recent introductions 
due to human activities (Othman et al. 2020) from a hith- 
erto unknown region of the Indian subcontinent (Fig. 2). 


Taxonomy 


All the hereby presented markers and overall phylogenetic 
pattern of the clade to which sequences of D. melanostic- 
tus from Pakistan belong to (see Figs 1—3) support previous 
studies discussing the unsettled taxonomy of this toad com- 
plex in its entire distribution (Wogan et al. 2016; Bisht et al. 
2021). The taxonomy of this group is difficult to resolve due 
to missing data from different regions, especially from the 
western part of the range, and similar morphology. In our 
study, we provide, for the first time, the evidence that the 
“hazarensis’ population from Pakistan belongs to the wide- 
ly distributed clade specified as Duttaphrynus sp. 1 (16S, 
sensu Bisht et al. 2021), or clade B (tRNAGIy-ND3; sensu 
Vences et al. 2017), respectively. It implies that Khan (2001) 
correctly observed morphological differences on these toads 
(see Introduction) but the lack of broader comparative ma- 
terial probably did not allow to compare specimens from 
Pakistan with other specimens, especially from the Indian 
subcontinent. The description, moreover, is formally inval- 
id, since the holotype lacks information about the name and 
location of the repository where it has been deposited (Art. 
16.4, ICZN 1999). Khan (2001) compared the Pakistani 
population with island populations from Taiwan (which, ge- 
netically, probably belong to Duttaphrynus sp. 2 clade in the 
16S tree, and clade G in the tRNAGIy-ND3 tree, respective- 
ly) and Sri Lanka from where genetic data are lacking. 
Following the recommendation of Dubois and Ohler 
(1999), 1.e. “Jf further studies demonstrate that a distinct 
species, close to Bufo melanostictus, occurs in Bengal, the 
name Bufo bengalensis Daudin, 1802 is available for it’, 
and based on the /6S phylogeny, the geographic distribution 
of Duttaphrynus sp. | clade (sensu Bisht et al. 2021), and 
the Principle of priority rule (ICZN 1999), we suggest that 
the name Duttaphrynus bengalensis (Daudin, 1802) comb. 
nov. (type locality “Bengale”, neotype MNHN 4967; for 
description see Dubois and Ohler 1999) is applicable for this 
clade. “Bengale” is equal to the area that is divided between 
Bangladesh and several states in neighbouring India (states 
West Bengal, Jharkhand, Bihar and the Cachar region of 
Assam matching the distribution of the Duttaphrynus sp. 1 
clade (Fig. 2A). According to our data, the clade is native 
to Bangladesh, India, Nepal and Pakistan, and apparently 
introduced to Indonesia, Malaysia, and Singapore. (see 
above and Fig. 2). Morphologically, the descriptions of 
D. m. hazarensis (Khan, 2001) and D. bengalensis (Daudin, 
1802) comb. nov. are almost identical in several aspects: 
the interorbital space is concave in both “hazarensis” and 
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D. bengalensis comb. nov., the space is larger than the 
width of the upper eyelid and internarial distance. Tympanic 
oblique, its diameter tn both these species 1s about two 
thirds of the eye diameter; the tympanic-eye distance 1s 
about less than one half of its distance. Moreover, the larger 
size of “hazarensis”’ (max. SVL of female up to 111.0 mm) 
and D. bengalensis comb. nov. (max. SVL of female 
neotype up to 133.2 mm) differentiates this species from 
D. melanosctictus (Schneider, 1799) specimens from India 
whose maximum SVL is 71.6 mm (lectotype ZMB 3462; 
Dubois and Ohler 1999). The snout length in D. bengalensis 
comb. nov. and “hazarensis’ is longer than the horizontal 
diameter of eye (vs. shorter in D. melanostictus). The canthus 
rostralis in D. bengalensis comb. nov. and “hazarensis” is 
sharp, with a ridge in comparison to D. melanostictus whose 
canthus rostralis is angular. The tympanum diameter in 
D. bengalensis comb. nov. and “hazarensis” is oval, oblique, 
its diameter about two-thirds of the eye diameter, whereas 
it is vertical and its diameter is more than half of the eye 
diameter in D. melanostictus. The tympanum-eye distance 
in D. bengalensis comb. nov. and “hazarensis” is about less 
than one half of its diameter, whereas it is less than one sixth 
of its diameter in D. melanostictus. Thus we recommend to 
assign the name D. bengalensis (Daudin, 1802) for the clade 
Duttaphrynus sp. 1/clade B for future studies. 
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Explanation note: Dataset of sequences used for the study. 
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